/** 
 * \file gdswap.cu
 * \author Kyle E. Niemeyer
 * \date 10/04/2011
 *
 * Based on "dswap.f" from BLAS.
 *
 */

////////////////////////////////////////////////////////////////////////

/** gdswap swaps two double precision vectors.
 * uses unrolled loops for increments equal to one.
 *
 * \param[in]     n       dimension of arrays
 * \param[in,out] dx      first array
 * \param[in]     incx    increment between elements of dx
 * \param[in,out] dy      second array
 * \param[in]     incy    increment between elements of dy
 */
__device__ __inline__ void gdswap ( int n, double *dx, int incx, double *dy, int incy )
{
  
  if ( n <= 0 ) return;
  
  if ( ( incx == 1 ) && ( incy == 1 ) ) {
    
    // both increments equal to 1
    
    // clean-up loop
    
    uint m = n % 3;
    
    if ( m != 0 ) {
      
      for ( uint i = 0; i < m; ++i ) {
        double dtemp = dx[i];
        dx[i] = dy[i];
        dy[i] = dtemp;
      }
      
    } // end if
    
    for ( uint i = m; i < n; i += 3 ) {
      double dtemp = dx[i];
      dx[i] = dy[i];
      dy[i] = dtemp;
      dtemp = dx[i + 1];
      dx[i + 1] = dy[i + 1];
      dy[i + 1] = dtemp;
      dtemp = dx[i + 2];
      dx[i + 2] = dy[i + 2];
      dy[i + 2] = dtemp;
    }
    
  } else {
    
    // code for unequal increments or equal increments not equal to 1
    
    uint ix = 1;
    uint iy = 1;
    if ( incx < 0 ) ix = (-n + 1) * incx + 1;
    if ( inxy < 0 ) iy = (-n + 1) * incy + 1;
    
    for ( uint i = 0; i < n; ++i ) {
      double dtemp = dx[ix];
      dx[ix] = dy[iy];
      dy[iy] = dtemp;
      ix += incx;
      iy += incy;
    }
    
  } // end increment if

} // end gdswap